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Limit Cycle Oscillations in Pacemaker Cells. 



II. Existing Models 
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Abstract — In recent decades, several mathematical models de- 
scribing the pacemaker activity of the rabbit sinoatrial node have 
been developed. We demonstrate that it is not possible to es- 
tablish the existence, uniqueness, and stability of a limit cycle 
oscillation in those models. Instead we observe an infinite num- 
ber of limit cycles. We then display numerical results from a 
new model, with a limit cycle that can be reached from many 
different initial conditions. 

Keywords — Sinoatrial node, electrical activity, heart, mathe- 
matical model, nonlinear dynamics. 

I. Introduction 

In elementary electrostatics it is well known that the relation 
between the voltage and the charge of a capacitor is 
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where v is voltage, C is capasitance and q is charge. Differenti- 
ating this equation with respect to time we obtain 



dq 
~dt 
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dv 
~dt 



(2) 



where 



dq 



is a current, and the sign of v is a matter of 



convention. In physiology this second relation (0) is being used 
to describe how the membrane potential (v) is changed when 
ions move across the cell membrane. Unfortunately we observe 
that this equation also is being used in some models where one 
in addition keeps track of the charge (q) concentrations inside 
and outside the cell Q. In those models voltage and charge are 
believed to be independent dynamic variables: first one deter- 
mines the voltage by integrating the membrane currents, then 
one determines the charge by integrating the same membrane 
currents. 

The purpose of this article is to point out that integrating the 
membrane currents once is enough. Voltage and charge cannot 
simultaneously be independent dynamical variables in a model, 
simply because of (Q). 

In order to visualize the drawbacks of treating voltage and 
charge as independent variables, we explore numerically the 
nonlinear dynamics of two different models describing the pace- 
maker activity of the rabbit sinoatrial node. The procedure is 
as follows: 

1. we integrate numerically the equations of motion for a suffi- 
ciently long time to detect a steady state, 

2. we change the initial conditions and repeat 1. 

First we display results from the Wilders et al. model a 
model that treats voltage and charge as independent variables. 
In that model it is thus possible to select an initial voltage and 
an initial charge independently. The dynamics of that model 
seems peculiar. An infinity of limit cycles is observed: each time 
we select new initial conditions a new limit cycle, corresponding 
to a new value of the constant of motion q — Cv, is found. This 
hampers the usefulness of the model. 

Second we display results from a new model of Endresen et al. 
(]] , where the voltage is not a dynamic variable. Here we cannot 
select an initial voltage independently of the initial charge, and 
only one limit cycle is observed. 
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The model of Wilders et al. J2j of the pacemaker activity of 
the rabbit sinoatrial node serves as an excellent example of the 
many models where the membrane potential is thought to be 
independent of the intracellular and extracellular charge concen- 
trations. In that model the equation of motion for the voltage 
is given by (|) 



dv 1 . . . . . . .... . 

"-77 = — -^(lb,Ca+«b,K+ib,Na + *Ca,L+JCa,T+lf+*K+lNa + «NaCa+*NaK) 

(3) 

There are fifteen dynamic variables in that model, the voltage v, 
the gating variables cZl, d^, ft,, /t, x, y, h, m, p, and the ionic 
concentrations [Ca]i, [Ca] re i, [Ca] up , [K]i, [Na]i. We want to 
determine how the long term dynamics in that model is changed 
when we change the initial conditions. To keep matters simple, 
we only change one initial condition: the initial intracellular 
concentration of potassium ([K]i); and we study the dynamics 
in two dimensions only: the phase space of v and [K]i. 
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Fig. 1 

Phase portrait of the Wilders et al. model. In (a) the trajectories of the 
oscillator and in (b) the limit cycle and the initial conditions (Pi) given by ( 
and [K]j = 140. 



Figure 1 displays the two-dimensional dynamics of the model 
with the initial conditions (dimensions skipped) 



v = -60.03 
d L = 0.0002914 
d T = 0.0021997 
/l = 0.9973118 
/t = 0.1175934 



x = 0.3294906 
y = 0.1135163 
h = 0.1608417 
m = 0.1025395 
v = 0.2844889 



[Ca] ; = 0.0000804 
[Ca] fe] = 0.6093 
[Ca] up = 3.7916 
[K] ; = variable 
[Na] ; = 7.5 , 



(4) 



and [K] ; = 140. In figure 1 (b) the point Pi denotes the v and 
[K]i coordinates of the initial conditions, and the closed loop 
at the bottom is the limit cycle. The trajectory of the model 
is displayed in 1 (a) where we observe that the model spirals 
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downwards from the point Pi to the limit cycle. 
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Phase portrait of the Wilders et al. model. In (a) the trajectories of the 
oscillator and in (b) the limit cycle and the initial conditions (Pi) given by (^J 
and [K] ; = 130. 



If the limit cycle in figure 1 (b) is unique it should be pos- 
sible to reach it from another initial condition. Let us try an 
initial condition below the limit cycle, and investigate whether 
the model spiral up towards it. We change the initial concen- 
tration of potassium from 140 to 130, leaving the fourteen other 
initial conditions unchanged. The result is displayed in figure 2. 
The model does not spiral upwards to the limit cycle in figure 
1, instead the model spiral downwards to a different limit cycle. 
We observed numerically a new limit cycle for each new initial 
value of [K]i, implying the existence of an infinite number of 
limit cycles. The model's fundamental flaw is clearly demon- 
strated. 

A New Model 

In a new model Q of the pacemaker activity of the rabbit sinoa- 
trial node, the membrane potential is determined by (111) 
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{[K]; - [K] + 2([Ca]i - [Ca] e ) + [Na]; - [Na] c } , (5) 



where q = FV {[K]; - [K] e + 2([Ca]s - [Ca] e ) + [Na]; - [Na] e } 
is the charge difference, F is Faraday's constant and V is cell 
volume. Here the ionic currents alter the concentrations which 
in turn alter the voltage, i.e. the physical quantities were calcu- 
lated in the following order: current i => charge q voltage v. 
The model has five dynamic variables, the gating variables x, h 
and the ionic concentrations [Ca]i, [K]i, [Na];, and we use the 
initial conditions: 



x = 0.9165 
h = 0.0000 
[K]j = variable 
[Ca]j = 0.004094141 
[Na] ; = 18.73322695 



(6) 



In this model we first notice that the initial value of [K]i, due 
to (js]), is not independent of the initial value of the voltage 
v. Thus changing [K]i changes v as is always the case when 
charging a capacitor (|ll). Second we notice that a tiny change 



in [K]i corresponds to a large change in voltage v, since the 
constant FV/C is large in most situations. 
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Fig. 3 

Phase portrait of the Endrcscn et al. model with the initial conditions (^}. In 

(a) the trajectory with [K]^ — 130.650 (Pi), in (b) the trajectory with 
[K], = 130.655 (P2), in (c) the trajectory with [K]; = 130.662 (P3), and in (d) 
the unique limit cycle and the three initial conditions Pi, P2 and P3 from (a), 

(b) and (c). 



In figure 3 we have displayed the simulation results from the 
model of Endresen et al. pj with three slightly (due to the 
large constant FV/C) different initial values of [K] ; : 130.650 
(a), 130.655 (b), and 130.662 (c). In figure 3 (d) the three initial 
conditions Pi, P2 and P3 all converge towards the same limit 
cycle. In an extensive numerical study we have not observed 
any physiological initial conditions that do not converge toward 
this limit cycle. In fact the same limit cycle can be reached 
when starting from the full equilibrium situation with equal 
intracellular and extracellular ionic concentrations Q. 

III. Discussion 

We have displayed numerical results from two types of mathe- 
matical models of the pacemaker activity of the rabbit sinoatrial 
node. The first type of model Q] showed an infinite number of 
limit cycles, the second type of model |j]] a limit cycle that could 
be reached from many different initial conditions. In order to 
avoid the drawback with an infinite number of limit cycles seen 
in the first type of models, we suggest that one should not treat 
membrane voltage (v) as a dynamic variable. Instead one should 
calculate the voltage using ||y), or at least select the initial con- 
ditions in agreement with 
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